﻿Imports DotSpatial.Topology, DotSpatial.Data, DotSpatial.Symbology, System.Threading.Tasks
#Const SingleThreaded = True

''' <summary>
''' Encapsulates routines to convert polygon or point shapefile to raster, determine IDW weights, and reassign grid values based on attribute values
''' </summary>
''' <remarks></remarks>
Public Class clsRaster
    Private lyrFeature As FeatureLayer = Nothing
    Private FieldName As String = ""
    Private RasterName As String = ""
    Private CellSize As Double = 50
    Private NumClosestPoints As Integer = 5
    Private MaxDistance As Double = 500
    Private lstNeighbors(,) As List(Of Tuple(Of Integer, Double)) 'indices of nearest points in feature layer and their weight (Nothing if a NoDataValue cell; 1 to NumClosestPoints otherwise)
    Private lyrRaster As Raster = Nothing

    Sub New(lyrFeature As FeatureLayer, FieldName As String, RasterFileName As String, CellSize As Double, NumClosestPoints As Integer, MaxDistance As Double)
        Me.lyrFeature = lyrFeature
        Me.FieldName = FieldName
        Me.RasterName = RasterFileName
        Me.CellSize = CellSize
        Me.NumClosestPoints = NumClosestPoints
        Me.MaxDistance = MaxDistance
    End Sub

    Public Function CreateRaster() As Raster
        lyrRaster = DotSpatial.Analysis.VectorToRaster.ToRaster(lyrFeature.DataSet, CellSize, FieldName, RasterName)
        ComputeWeights()
        UpdateGrid()
        Return lyrRaster
    End Function

    Private Sub ComputeWeights()

        If lyrFeature.DataSet.NumRows = 0 Then Exit Sub

        With lyrRaster

            ReDim lstNeighbors(.NumColumns - 1, .NumRows - 1)

#If SingleThreaded Then
            For r As Integer = 0 To .NumRows - 1
#Else
            Parallel.For(0, .NumRows, Sub(r As Integer)
#End If

                For c As Integer = 0 To .NumColumns - 1
                    If .Value(r, c) <> .NoDataValue Then
                        Dim x As Double = .Xllcenter + c * .CellWidth
                        Dim y As Double = .Yllcenter + r * .CellHeight
                        Dim pt As New Coordinate(x, y)
                        Dim dictDist As New SortedDictionary(Of Double, Integer)
                        For f = 0 To lyrFeature.DataSet.Features.Count - 1
                            Dim pc As Coordinate = lyrFeature.DataSet.Features(f).Centroid.Coordinates(0)
                            Dim dist As Double = pt.Distance(pc)
                            If dictDist.Count = 0 Or dist < MaxDistance Then dictDist.Add(pt.Distance(lyrFeature.DataSet.Features(f).Centroid.Coordinates(0)), f)
                        Next

                        lstNeighbors(c, r) = New List(Of Tuple(Of Integer, Double))

                        Dim sumwt As Double = 0.0

                        For Each kv In dictDist
                            If lstNeighbors(c, r).Count = 0 Or kv.Key < MaxDistance Then
                                Dim wt As Double = 1.0 / kv.Key 'inverse distance
                                sumwt += wt
                                lstNeighbors(c, r).Add(Tuple.Create(kv.Value, wt)) 'add index and weight if closest one or within radius of influence
                                If lstNeighbors(c, r).Count = NumClosestPoints Then Exit For
                            End If
                        Next

                        'recalculate the weights
                        If sumwt > 0 Then
                            For Each kv In lstNeighbors(c, r)
                                kv = tuple.Create(kv.Item1, kv.Item2 / sumwt)
                            Next
                        End If
                    End If
                Next

#If SingleThreaded Then
            Next r
#Else
                                      End Sub)
#End If
        End With
    End Sub

    Public Sub UpdateGrid()

        If lyrFeature.DataSet.NumRows = 0 Then Exit Sub

        With lyrRaster

#If SingleThreaded Then
            For r As Integer = 0 To .NumRows - 1
#Else
            Parallel.For(0, .NumRows, Sub(r As Integer)
#End If

                For c As Integer = 0 To .NumColumns - 1
                    Dim v As Double = 0
                    If lstNeighbors(c, r) Is Nothing Then
                        v = .NoDataValue
                    Else
                        For Each t In lstNeighbors(c, r) 'item 1 is index, 2 is weight
                            Dim f As Feature = lyrFeature.DataSet.Features(t.Item1)
                            v += TestNull(f.DataRow.Item(FieldName), 0) * t.Item2
                        Next
                    End If
                    .Value(r, c) = v
                Next

#If SingleThreaded Then
            Next r
#Else
                                      End Sub)
#End If
        End With
    End Sub
End Class
